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Abstract : 

Numerical results are presented for a generalized resonance valence bond state 
which includes both ionic and covalent contributions in each bond ; non nearest-neighbor 
bonds are also considered. Variational calculations have been performed and the space 
group symmetry has been taken into account. The results for the Hubbard Hamiltonian 
are compared with the exact ones for the 2x2x2 and the 4x4 lattices. The agreement 
is quite satisfactory for large values of the interaction U ( for U=40 the overlap with 
the exact ground state wave function is 0.999 ). Several correlation functions are also 
compared. 

PACS numbers: 75.10.Jm, 740.Mg 



1. Introduction. 



It is still under debate if the 2D one-band Hubbard model can provide the basic 
physics of high temperature superconductivity. Independently of this, and in spite of 
the obvious oversimplification introduced by the model in describing a strongly corre- 
lated electron system, the theoretical challenge offered by the 2D Hubbard model is 
certainly full of fascinating and unsolved problems. Only rough and tentative phase 
diagrams concerning antiferromagnetism, ferromagnetism and superconductivity have 
been produced, manly in the large U case, where the Hubbard model maps into the 
t — J model. 

It is interesting to compare the situation of our present knowledge of the Hubbard 
model and the theory of the fractional quantum Hall effect (FQHE). Both theories 
refer to strongly interacting two-dimensional electron systems; however the situation is 
much more satisfactory in the FQHE case, not only because of the simplicity of the 
experimental results, but for another fundamental reason: in the FQHE case we have at 
our disposal a valid theoretical tool, provided by the Laughlin's wave function [1]. All 
subsequent theories, like Jain's theory [2], make use of this beautiful form of the wave 
function. 

In the case of the Hubbard model, we do not have a form of the wave function quite 
as satisfactory. 

Many years ago, Gutzwiller [3] proposed a wave function obtained by applying a 
suitable operator to the Fermi sea wave function. The physical idea is simple: one tries 
to reproduce the main effect of the interaction, by weighting with a factor i] D all the 
configurations containing D doubly occupied sites; r/ is a variational parameter. The 
interesting rj = case corresponds to a state with no doubly occupied sites (Gutzwiller 
projection). 

The Gutzwiller wavefunction ipo provides a noticeable insight into the physics of 
the Hubbard model. However, the expectation values with respect to ipo cannot be 
computed without making some approximations. Furthermore the value of the ground 
state energy, tested in finite size cases, is not more accurate than the one given by the 
Hartree-Fock solution (spin density wave). 

In 1987 Anderson [4] proposed that the electronic wave function of a high T c su- 
perconductor is a "resonating valence-bond" (RVB) state. This state is obtained as the 
Gutzwiller projection of a BCS state with definite particle number. It can also be con- 
sidered as a superposition of valence bond states, i.e. states describing pairs of electrons 
in different sites coupled to a spin singlet. A "covering" of the lattice is a state in which 
all the electrons are coupled in pairs. 

A number of calculations (see, e.g. [5] ) have been performed on the RVB state, 
especially in a simplified version in which only singlets formed by nearest-neighbor 
electrons are allowed (NNRVB). Almost all these works refer to the Heisenberg and t-j 
two-dimensional models. Calculations involving NNRVB states are made possible by 
simple graphic rules [6] , [7] . 
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In this type of approach the Gutzwiller projection is not necessary, since the 
coverings of a lattice at half-filling are products of dimers created by operators like 
Cij = c t\ c tl c jT c iI > wnere * an d 3 are nearest- neighbor sites and c^,c^ denote the 
usual creation operators; therefore no doubly occupied site occurs. 

It is interesting to notice that long before Anderson's proposal of the RVB state 
for strongly interacting electrons, some theoretical chemists, Choi and Thorson [8] com- 
puted the one and two-particle density matrices for a generalized RVB state, both for 
a linear chain and a 4 x 4 lattice. By " generalized RVB state " we mean an RVB 
state such that the operator + creating the valence bond between the sites i and j 
contains not only a covalent term C^, but also a ionic contribution 1^ = c^c^ + c ^- c tp 
therefore 0+ = ACj + where A is a variational parameter. These authors 

mentioned the Mott insulator and the Hubbard Hamiltonian as possible applications, 
but no explicit calculations of the ground state energy was presented. 

This generalized RVB state has been proposed again for the Hubbard model, and 
independently, by Messager and Richard [9] ; it was also studied for a strip of 2 x N sites 
by Cicuta and Stramaglia [10]. To the best of our knowledge, there has been no further 
study on this proposal, despite the following arguments in favor of its consideration as 
a good variational state for the Hubbard model: 

1) The state O + |0 > is the exact ground state for a two electron system. 

2) The state (O^O^ — 0^0^)10 > is an extremely accurate approximation of the 
exact ground state of a 2 x 2 Hubbard model. ( We recall, that in the analogous case 
of the FQHE, the analytic part of the Laughlin's wave function is a product of factors 
(zi — Zj) m , which correspond to the exact ground state for two particles; therefore the 
idea of building a many-particle wave function by combining some kind of product of 
two-particle wave functions deserves particular attention. Notice that the three particle 
Laughlin's w.f. is exact on the sphere geometry.) 

3) The introduction of a parameter A that controls the amount of double occupancy 
of the bond seems to be a natural way of constructing a wave function which covers 
both the large U case ( A — > 1 ) and the small or medium U case. 

4) The ground state energy of the large U limit of the Hubbard model at half filling 
( i.e. the Heisenberg model ) can be reproduced with great accuracy by including in 
the RVB state also bonds connecting sites which are far apart. Indeed Liang, Doucot 
and Anderson [11] have studied this state on a bipartite lattice, by allowing bonds 
connecting only sites i of the sublattice A with sites j of the sublattice B, and giving a 
weight factor h(\i — j\) to each bond, where the factor h decreases with a power law with 
the distance \i — j\; notice that the positivity of h guarantees the important Marshall 
sign rule [12] . 

5) Even if we do not succed in finding a form of the wave function as useful as 
Laughlin's w.f., we believe that a good variational w.f. can be of use as an initial guess 
for modern numerical methods. 

Numerical calculations on the Hubbard model can be performed either with great 
accuracy on small systems, or using Montecarlo or other approximations on larger sys- 
tems. Here we will be concerned with exact results on small systems (2x2x2 and 
4x4). Full configuration interaction calculations at half filling on a 4 x 4 lattice are 
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greatly simplified by the consideration of the space symmetry group of the model. The 
best way of taking into account the space symmetry group is to replace the electron 
configurations by linear combinations of configurations contained in a group orbit ( we 
recall that a group orbit is the set {g\c >,g G G} where \c > denotes a configuration 
) ; the coefficients of the linear combination are determined by the irreducible repre- 
sentations of the group, and are constant for the totally symmetric state. The orbits 
are of variable lenght; for instance the Neel configuration \n > gives rise to a small 
orbit, constituted by \n > itself and the "anti-Neel" configuration \n > with reversed 
spins ( indeed all rotations and translations of the lattice either send \n > — > \n > or 
| rc > — > \n > ). The dimension of the Hilbert space can be greatly reduced since the 
order \G\ of the space group is large . For a square n x n lattice with periodic boundary 
conditions the translation group in the two directions has n 2 elements, reflections in 
the direction of the main diagonals and one of the axis produce a further factor of 8; 
this gives \G\ = 128 for n=4 ( Actually the 4x4 lattice is topologically equivalent to 
an ipercube in 4 dimensions, so a further symmetry arises and a larger group G can be 
found such that \G\ = 384 ). For a group-theretical analysis of the Hubbard model see 
ref.[13]. 

Besides space symmetry, the Hubbard Hamiltonian presents the usual spin symme- 
try. It is a very difficult task to use an orthogonal basis adapted to both symmetries. 
One of the advantages of the RVB approach, is that both space and spin symmetries 
are introduced in the wave function from the outset. Unfortunately, like in all valence 
bond methods, the wave functions to be used are not orthogonal. 

In the present work we have considered 2x2x2 and 4x4 systems with periodic 
boundary conditions and we have computed the expectation value of the energy and 
other observables (spin-spin and pairing correlation functions ) for a generalized RVB 
state. Since on a bipartite lattice the next type of bond to be considered beside the 
nearest- neighbor case is provided by bonds between electrons at distance \/h on the 
lattice ( named KM or chess knight's move bond [14] ), we have taken into considerations 
wave functions admitting also one or more KM bonds. 

All the numerical results presented here correspond to half-filling. However, the 
formalism can be easily extended in order to cover differnt fillings (see the end of Sec. 
2). 

2. Formalism . Let us consider an M x N square lattice. A site i is determined by 
two integers (r, s). The sublattice A is formed by the sites such that (— l) r+s = 1, and 
the sublattice B is formed by the sites such thtat (— l) r+s = — 1. Only bonds connecting 
the sublattices A and B will be allowed. 

We define now the creation operator of a generalized valence bond between the 
sites i and j: 

(1) A + (x,y) = xl+ + yC+ 

where x and y are variational parameters and if^C^ are defined by: 

( 2 ) % = ctA + c fcti 



3 



(3) C ± - c+ c+ + c+ c+ . 

The operator C% creates the usual covalent bond, while 1% creates an ionic bond. 
Normalization of the state vector Ajj(x,y)\0 > implies the condition: 

(4) 2x 2 + 2y 2 = l 

A covering of the lattice is obtained by considering the following product of oper- 
ators applied to the vacuum: 

(5) \(f> k >=A+ ji (x u y 1 )A+ h (x 2 ,y2)...A+ jn (x n ,y n ) |0 > 

where i±,Z2, ■ ■ - in belong to the sublattice A , ji, j 2 , . . .j n belong to the sublattice B, 
and the index k runs over all possible coverings. 

At half-filling the number of lattice sites is twice the number of dimers, i.e. NM = 
2n. Because of condition (4), the states \cf>k > are normalized. However they are not 
orthogonal, and it is also possible that they are linearly dependent. 

We assume that in each factor A^Ax^y) the variational parameter x is a function 
of the distance \i — j\. 

A generalized RVB state is a linear combination of coverings: 

(6) H>= J2 Ck ^ k > 

k 

Let g denote an element of the space group. We assume that the coefficients cj~ are the 
same for all the coverings g\(f>k >, g £ G. In this way \ip > becomes a totally symmetric 
state. It must be noticed that usually the ground state of a finite size Hubbard model 
is totally symmetric with respect to the space group; an exception is provided by the 
2x2 Hubbard model. 

In order to compute the matrix elements of the overlap matrix Sij =< 4>i\<pj >, 
Sutherland's graphical rules [6] can be used. The scalar product < <f>i\<j)j > can be 
written as < 4>i\4>j >= YlT=i ^ where l k denotes the contribution of a loop with an 
even number of bonds. In order to compute the factors the following simple algebraic 
rules can be employed. From the commutator: 

(7) [A ij (x 2 ,y2),A±(x 1 ,y 1 )] = 2(x 2 x 1 + y 2 y 1 ) 
-{x 2 x 1 + y 2 yi)(c+c^ + cjc a + c+ c jT + c+ c n ) 

-{X2V1 + V2 x l)( c t^ + c t c jl + S + t Ci T + C+jC^) 
we can deduce the formula: 

< X\A ij (x 2 ,y 2 )A+ j (x 1 ,y 1 )\Y >= 2(x 2 x 1 + y 2 y 1 ) < X\Y > 

for any partial coverings \X >,\Y > such that the sites i,j are not occupied (i.e. 
coverings of the lattice without the sites i and j ) . 
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Analogously, from the commutator (k ^ j): 

(9) [A ik (x 2 ,y 2 ) 1 Af j (x 1 ,y 1 )] = x 2 x x -x 2 xi(c+ t q t +cjc a ) -x 2 yi(c^c^ + c^a^ 

-^2^i(c, + T CfcT + 4l c H) ~ ^22/i(S + T Cfc T + c tl c kl) 

we obtain: 

(10) < X\A j i(x 3 ,y 3 )A ik (x 2 ,y 2 )A±(x 1 ,y 1 )\Y > = < X\Ai k (x!X 2 xs, -y^ysW > 

where in the partial coverings \X >, \ Y > the sites i,j, k, I are empty. 

In Fig. 1 it is shown how formulas (8), (10), can be used to simplify a Sutherland's 
graph. 

From relations (8), (9), it follows that the contribution of a loop with m nearest- 
neighbor bonds is given by (we assume x, y real): 

(11) l = 2x m + (-l) 1+ ^2y m 

in agreement with previous results [10], [11]. 

Let us now denote by A^-(x,y) the operator creating a "knight's move bond" ( 

KM bond) between the sites i and j whose relative distance is y/b (l 2 + 2 2 = 5)!. The 
introduction of a KM bond simply modifies the parameters x, y into x,y; therefore the 
contribution of a loop containing r KM bonds and m — r nearest-neighbor bonds is given 
by: 

(12) I = 2x m ~ r x r + {-l) 1+ ^2y m - r y r 

Formulas (11), (12) allow us to compute the overlap matrix < <pi\cj)j >. 

Let us now compute the matrix elements < (f) r \H(f> s > of the Hubbard Hamiltonian. 
Denoting by the operator: 

(13) Tij = c+c jT + cjcjj. + c+ c lT + c+ c a 
the Hubbard Hamiltonian can be written as: 

(14) H = - J2 T u ■ ''»'••»'•. 

<ij> i 

where '^2 <i j > denotes the sum over all nearest-neighbor pairs. 

As for the matrix elements of the operator T^, we must distinguish three cases: 

i) The operator TV, acts on a pair of sites i,j belonging to two different closed 
loops. Since the total particle number must be conserved inside a loop, the contribution 
vanishes. The corresponding graph is shown in Fig. 2 . 

ii) The operator acts on a pair of sites connected by a bond. For any partial 
coverings \X >, \ Y > such that the sites i, j are empty, we have : 

(15) < X\T t3 At 3 (x,y)\Y >= 2 < X\A±(y, x)\Y > 
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Hence the hopping operator exchanges the covalent and ionic bonds. 

iii) The operator T^j divides a loop into two parts of lengths m, n; if m and n are 
odd, by repeated application of Eq. (10) we reduce ourselves to the preceding case, and 
Eq.(15) still holds; if m and n are even we can repeat the same procedure. In this way 
we arrive to an expression of the form: 

(16) < 0\A il (x 1 ,y 1 )A jm (x2,y2)T lJ A+(x3,y3)A+ n (x 4 ,y 4 )\0 > 

whose contribution vanishes. The corresponding graph is shown in Fig. 3. 

Let us denote by Vij = U(n^riii + rij^nji) the contribution to the potential energy 
of the sites % and j. Since Vij counts the double occupations, the action of Vij on 
Af- (x, y) is simply to annihilate the covalent part of the bond. Hence: 

(17) < X\V tJ A+(x,y)\Y >= U < X\A+(x,0)\Y > 

Let us add a few words on the modifications of the formalism that are necessary to 
cover cases different from half-filling. In general, the scalar product of two coverings 
can still be easily computed. However, the formula < 4>i\4>j >= YlT=i ^ does not hold 
any more. The Sutherland's graph corresponding to < <f>i\<f>j > contains in general not 
only closed loops, but also open paths. These paths can be easily simplified by using 
repeatedly formula (10). If the number of bonds of a given open path is odd, the final 
result vanishes, since < Aij(x,y) >= ; if it is even, the final contribution is of the 
type : 

< A ij j(xi,yi)Af k (x2,y2) >=xlx 2 . 

Notice that only the ionic coefficients x\,X2, enter in the last formula . 

The generalization of the rules to be used in order to compute < <f>i\H\<f)j > is 
straightforward. 

Numerical results and conclusions. 

In order to obtain a first insight about the values of the variational parameters, the 
scalar product with the exact ground state, etc., we have first considered a toy model, 
i.e. the 2x2x2 cube, which is topologically equivalent to the 4x2 lattice (see e.g. ref. 
[14] for the case of the Heisenberg model ). 

Denoting the vertices of the cube by the numbers 1,2, .. .,8, we can choose five 
coverings that give rise to five different orbits of the space group; in Fig. 4 the coverings 
that generate these five orbits are shown: there are two coverings without KM bonds, 
and three coverings with 1, 2 and 3 KM bonds respectively. 

Since we are interested in the totally symmetric state, all the covering structures 
of an orbit are taken with equal weight. 

Notice that the KM correspond to bonds between opposite vertices of the cube. 
The matrix elements S rs , T rs , V rs (r, s = 1, 2, ..5) of the overlap, hopping and interaction 
matrices can be easily computed (see Appendix) , where the indexes r, s run over the 
five orbits of Fig. 4 . 
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In Table I we give, for the 2x2x2 case, the value of the energies of the knight's 
move resonating valence bond state |KM > and the nearest neighbor resonating valence 
bond state |NN >, compared with the exact energy E and the best Hartree-Fock energy 
Ehf; we denote by |E > the exact ground state and we give the values of the overlaps 
< E|NN >, < E|KM >. 

Even in this toy model we can learn something: the RVB variational states (nearest- 
neighbor and knight's move) are worse than the Hartree-Fock state for small U, but are 
superior (especially the KMRVB state) for U > 8. 

We have repeated the same calculation for the more interesting case of the 4x4 
lattice. For this model there are 272 NNRVB coverings; this number increases to 40320 
if we include also KMRVB coverings. In Fig. 5 an exemple of covering of the lattice 
containing 3 KM bonds is shown. 

We know from previous works (see, e.g. ref.[13]) that the ground state is totally 
symmetric with respect to the space group G. Averaging over G, we can reduce the 
dimension of the secular problem to 13 in the nearest-neighbor case, and to 458 if 
KM bonds are present. For comparison, the dimension of the whole Hilbert space is 
(g) ~ 1.6 10 8 , and the dimension of the totally symmetric subspace is ~ 1.3 10 6 . 

In Table II the values of the ground state energies and the overlaps with the ground 
state are given for the 4x4 case. The notations are the same as those used in Table I. 

From Table II we see that for large U (e.g. for U=40) the KMRVB state has an 
energy value (-1.8973 compared to -1.9084) and an overlap with the exact ground state 
of 0.999. Both are very satisfactory values. 

It is also interesting to remark that considering the projection of the exact ground 
state wave function \E > in the subspace generated by the coverings \4>k >, the absolute 
values of the coefficients =< <j>k\E > varies very little ( typical values are 0.025, 
0.0248, 0.0253, 0.0269 ); furthermore their sign agrees with the Marshall-Peierls rule; 
this kind of behaviour is well known for the Heisenberg model (see ref. [14]). The 
variational state considered in refs.[10],[ll] consists in taking all the \ck\ equal; this 
approximation is very good in the restricted class of the NNRVB states. 

In Table XXXX the variations of the ionicity parameters r = - , f = S with U are 



exibited for the models 2x2x2 and 4x4. 

From the Table we see that the values of the variational parameters are not exactly 
the same for a given U, but their order of magnitude is the same. Hence we can speak 
of an " approximate transferability " from one model to another. 

In table III we show the results for the spin-spin correlation function < S z (0)S z (r) > 
for r = (0, 0), (1, 0), (2, 0), (2, 1), (2, 2); (in order to save space, in the table only the value 



As is well known the antiferromagnetic order of the Hartree-Fock state is too strong 
compared with the exact ground state values, since the best mean field solution is a spin- 
density wave; the antiferromagnetic order of the NNRVB state is reasonable but decays 
too fast as expected (see also ref [8] ). The spin-spin correlation function derived from 
the KMRVB state is almost exact for U=40, but even for U=4 it reproduces the correct 
behaviour, while this is not the case for the Hartree-Fock solution. 

In table IV we have examined the d-wave pair correlation function in order to see 
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of |r| = \/r% + r% is indicated). 
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the effect of the approximation on a physical quantity which is probably crucial in high 
T c superconductivity. This correlation function is defined as < A+(0)A(r) > where: 



and denoting by x and y unit vectors in the two lattice directions, /(±x) = 1 and 
f(±y) = — 1. Both correlation functions, derived from the NNRVB and KMRVB wave 
functions, show a quite reasonable behaviour for the rate of decay at increasing values of 
|r|. We have also computed the s-wave pair correlation function which shows a similar 
behaviour. This kind of agreement for the s and d-wave pair correlation functions is 
present only for large values of U. No indications of ODLRO appears. 

Concluding, the type of variational wave functions studied in this work certainly 
captures a good part of the physics of the Hubbard model, especially for large values of 
U. Hence, it would be of interest to study the excited states, obtained by substituting 
one triplet bond for a singlet, or to introduce holes in the system. A strong limitation 
of the formalism is constituted by the non orthogonality of the coverings. 





u=±x,±y 
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Appendix 

In the following, the matrix elements S rs , T rs , V rs (r,s = 1, 2, ..5) of the overlap, 
hopping and interaction matrices for the 2x2x2 cube are given, (x, y) are the parameters 
of KM bonds, (x, y) are the parameters of NN bonds. 

OVERLAP MATRIX 

5 11 = \[l + 2(x 2 -y 2 ) 2 ] 

512 = \{x 2 - v 2 )[^ + y 4 + 2] 

Sis = l[6(x 5 x + y 5 y)] 

S 1 ,= 1 -[2(x 2 x 2 -y 2 y 2 )+4(x 6 x 2 -y 6 y 2 )] 

S 1B = [2(x 2 x 2 - y 2 f)f 

522 = \[l + (x 2 - y 2 ) 2 + (x 6 + y 6 )] 

523 = ^[Q(x s x~y 3 y] + 6(x 7 x - y 7 y)] 

524 = l[S(x 4 x 2 + y 4 f) + 8(x 3 x - y 3 y) 2 \ 

5 25 = 2(x¥- !/ V) 

533 = i [1 + 6(x 4 x 2 + y 4 f) + 2(x 6 + y 6 ) + 12(x 3 x - y 3 yf] 

5 34 = (x 3 x - y 3 y) + (x 5 x 3 - y 5 y 3 ) 

5 35 = 2(x 3 x 3 + y 3 y 3 ) 

Sm = \[1+ [2(x 2 x 2 - y 2 y 2 )] 2 + 8(x 4 x 2 + y 4 y 2 )] 
S 45 = 2(x 2 x 2 -y 2 y 2 ) 
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$55 — 1 

HOPPING MATRIX 

T 11 = ^xy(l + 2(x 2 -y 2 ) 2 ) 

T 12 = -2- xy[2(x 6 - y 6 ) + (xy) 2 (y 2 - x 2 ) + 4(x 2 - y 2 )] 

T i3 = ^[16xy(x 5 x + y b y) + 20xy(x 3 x + y 3 y) - 8x 2 y 2 (yx + xy)] 

T14 = ^{2[24xy(x 4 x 2 - y 4 f) + 16xyx 2 y 2 (y 2 - x 2 ) + 8x 3 y 3 (y 2 - tx 2 )] + 
+ 32xy(x 2 x 2 - y 2 f) + 8xy(x 2 - y 2 )} 

T 15 = 32[xy(x 2 -y 2 )(x 2 x 2 -y 2 y 2 )] 
1 fi 

T 22 = y [1 + 4(x 4 + y 4 ) - 3(xy) 2 + 2(x 6 + y G )] 

T23 = ^[32xy(x s x - y 3 y) + 12xy(xx - yy) + 

+ 28xy(x 5 x — y 5 y) + 12x 2 y 2 (—x 3 y + y 3 x) + 8x 3 y 3 (—xx + yy)] 

T 24 = - [2{lbxy{x 4 x 2 - y 4 f) + 16xy(x 2 x 2 + y 2 f) - 16x 2 y 2 xy+ 
+ xy(x 2 x 2 + y 2 y 2 )] + 48xy(xx - yy)(x 3 x - y 3 y)] 

T 25 = 16[xy(x 2 x 4 - y 2 y 4 ) + xyx 2 y 2 (y 2 - x 2 ) + x 2 y 2 xy(y 2 - x 2 )] 

T 33 = - [24xy + 2Axy{x 4 + y 4 ) + 3[16xy(x 4 x 2 + y 4 y 2 ) + 16xy(x 2 x 2 + y 2 y 2 ) + 
8 

+ 4x 2 y 2 xy] + 3 • 48xy(xx — yy)(x 3 x — y 3 y)] 
1 

T34 = - [12xy(xx — yy) + 16xy(x 3 x — y 3 y) + 20xy(x 3 x 3 — y 3 y 3 ) + 
+ 4xyxy(-x 3 y - y 3 x) + 24x 2 y 2 xy(-xx + yy)] 

T 35 = 12[xy(xx 3 + -yy 3 ) - xxyy(xy + yx)] 
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T 44 = - [I6xy + 4(16xy(x 2 x 2 + y 2 y 2 ) - \Qx 2 y 2 xy) + 32xy(x 2 x 2 - y 2 y 2 )(x 2 - y 2 )} 
T 45 = 8xy(x 2 - f) 
T55 = 

INTERACTION MATRIX 

V 11 = ^[l + 2x 2 (x 2 -y 2 )] 

V l2 = ^x 4 (l + 5x 4 ) 

y 13 = 4x 2 (x 5 x + y 5 y) + 6x 5 oc 

V 14 = -[16x 6 x 2 + 8x 2 (x 2 x 2 - y 2 y 2 ) + 4x 2 x 2 } 
V 15 = 16x 2 x 2 (x 2 x 2 -y 2 f) 
V 22 = ^[1 + 2(x 6 + y 6 ) + 6x 4 + x 2 (x 2 - y 2 )] 
V 23 = \[8x 2 {x 2 x 2 - y 2 y 2 ) + 4x 2 x 2 + 8x 7 x] 

F 24 = hl2x 4 x 2 + 8x 2 {x 4 x 2 + y 4 f) + 16x 3 x(x 3 x - y 3 y)} 
V 25 = 16x 4 x 4 

V 33 = \[12x 2 {x^x 2 + y 4 y 2 ) + 9x 4 x 2 + 48x 3 x(x 3 x - y 3 y) + 
8 

+ 6x 6 + 4x 2 (x 6 + y 6 ) + 6x 2 + 2x 2 ] 
V34 = 2x 3 x + 2(x 2 + x 2 )(x 3 x - y 3 y) + 8x 5 x 3 
V 35 = 4x 2 (x 3 x 3 + y 3 y 3 ) + 6x 3 x 3 

V M = -[x 2 + 2x 2 + 12x 4 x 2 + 8x 2 (x 4 x 2 + y 4 f) + 8x 2 x 2 (x 2 x 2 - y 2 y 2 )] 
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V 45 = 8x 2 (x 2 x 2 - y 2 f) + Ax 2 
V 55 = 8x 2 
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Table I. Energy values of the ground state, Hartree-Fock and RVB states for the 
2x2x2 lattice. 



u 


E 


E HF 


En n 


Ei KM 


(E\NN) 


(E\KM) 





-12.0000 


-12.0000 


— 8.8000 


-9.3463 


0.5477 


0.6941 


1 


-10.1187 


-10.0000 


— 7.2887 


-7.7171 


0.5097 


0.6765 


4 


-5.9542 


-5.0274 


-4.8152 


-5.0079 


0.6887 


0.8101 


8 


-3.4671 


-2.8341 


-3.1667 


-3.2697 


0.9218 


0.9763 


16 


-1.8772 


-1.4772 


-1.8002 


-1.8464 


0.9760 


0.9986 


40 


-0.7763 


-0.5985 


-0.7553 


-0.7728 


0.9850 


0.9993 



Table II. Energy values of the ground state, Hartree-Fock and RVB states for the 
4x4 lattice. 



u 


E 


Ehf 


Enn 


Erm 


(E\NN) 


(E\KM) 


1 


-20.7936 


— 20.6542 


-15.5266 


-16.9183 


0.44317 


0.64978 


4 


-13.6218 


— 12.5665 


-10.7020 


-11.5642 


0.58962 


0.79090 


8 


-8.4689 


-7.3896 


-7.2782 


-7.7707 


0.78905 


0.94436 


16 


-4.6119 


-3.9116 


-4.2421 


-4.4830 


0.89713 


0.99489 


40 


-1.9084 


-1.5941 


-1.8051 


-1.8973 


0.92734 


0.99904 
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Table III. Values of the spin-spin correlation function for the exact ground state, 
the Hartree-Fock and the RVB states. 



(Sz(0)S s (r)) 




r 


HF 


NNRVB 


KMRVB 


GROUND 


U = 4 





0.18704 


0.21316 


0.20883 


0.19244 


1 


-0.13706 


-0.08782 


-0.08841 


-0.06879 


2 


0.12321 


0.03401 


0.05120 


0.03750 


V5 


-0.12501 


-0.02010 


-0.05170 


-0.04520 




0.12201 


0.01447 


0.04442 


0.03851 


U = 40 





0.24876 


0.24863 


0.24855 


0.24853 


1 


-0.24784 


-0.11094 


-0.11579 


-0.11590 


2 


0.24753 


0.05014 


0.07110 


0.07085 


V5 


-0.24753 


-0.03262 


-0.06811 


-0.06750 


2^2 


0.24753 


0.02476 


0.06048 


0.05997 



Table IV. Values of the pair correlation function for the exact ground state, and 
the RVB states. 



(At(O)A(r)) d-wave, U = 40 


r 


NNRVB 


KMRVB 


GROUND 





2.33389 


2.39191 


2.38947 


1 


0.59034 


0.60339 


0.60118 


2 


0.00199 


0.00114 


0.00252 


V5 


0.00002 


0.00003 


0.00005 


2V2 


0.00000 


0.00000 


0.00000 


V2 


0.00589 


0.00487 


0.00261 
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FIGURE CAPTIONS 
Fig.l The operators Aij (dotted line)and Af- (full line) are represented by a bond con- 
necting the points i and j. 

Fig. 2 Graph representing the operator TV, acting on a pair of sites i,j belonging to two 
different closed loops. The corresponding contribution vanishes. 

Fig. 3 Graph corresponding to the expectation value : 

< 0\A il (xi,y 1 )A jrn (x2,y2)T ij A+(x 3 ,y 3 )Af m (x4,y4)\0 > 
whose contribution vanishes . 



Fig. 4 Representative coverings of the 2x2x2 lattice which generate different orbits 
of the space group. 

Fig. 5 An exemple of covering of the 4x4 lattice containing 3 KM bonds; point 13 and 
point 16 are connected by a bond because of periodicity. 
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Fig. 4 
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